################################################################################
#                                                                              #
#        Replication script for analysis contained in:                         #
#        Does the Presidency Moderate the President?                           #
#        Presidential Studies Quarterly (Vol. 47, Issue 1, 2017).              #
#                                                                              #
#        Barry C. Edwards, J.D., Ph.D.                                         #
#                                                                              #
################################################################################

########## Set-up for analysis ###########################
#  Remove all objects just to be safe
rm(list=ls(all=TRUE))
# library(foreign)
# library(gdata)

#### Change the working directory in line below for your own machine.
setwd('C:/Users/Barry Edwards/Dropbox/does presidency moderate president')

presModerationData <- read.csv("presidents who served in congress.csv", fileEncoding = "latin1")
# names(presModerationData)
attach(presModerationData)
# detach(presModerationData)


########################### Figure 1 summarizes hypotheses

# Uncomment next line and the corresponding dev.off() call to generate .png file images
# png(file="moderationHypotheses.png",width=4,height=4, units="in", pointsize=10, res=600)

par(mar=c(3.1,3.1,0.6,0.6), omi=c(0.1,0.1,0.1,0.1), family="serif", mfrow=c(1,1))
plotSize <- 1
plot(x="",y="", xlim=c(-plotSize,plotSize), ylim=c(-plotSize,plotSize),
     ylab="",xlab="", axes=F)
mtext(side=1,text="<< More Liberal                          More Conservative >>",cex=.9,line=.5,adj=.5)
mtext(side=1,text="Preferences as Member of Congress",cex=1,line=2,adj=.5,font=2)
mtext(side=2,text="Preferences as President",cex=1,line=2,adj=.5,font=2)
mtext(side=2,text="<< More Liberal                          More Conservative >>",cex=.9,line=.5,adj=.5)

polygon(x=c(0,2,2), y=c(0,0,2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(0,-2,-2), y=c(0,0,-2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(0,0,2), y=c(0,2,2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
polygon(x=c(0,0,-2), y=c(0,-2,-2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
abline(h=0,lty=1)
abline(v=0,lty=1)
abline(a=0,b=1,lty=1,lwd=8,col="gray20")
text(x=c(.7,-.7), y=c(.3,-.3), label="Moderated\nPreferences", cex=.8)
text(x=c(.3,-.3), y=c(.7,-.7), label="Amplified\nPreferences", cex=.8)
text(x=c(.55,-.55), y=c(-.55,.55), label="Not Expected", cex=.8)
text(x=c(.55,-.55), y=c(.55,-.55), label="NO IDEOLOGICAL CHANGE",srt=45, cex=.6, col="white",font=2)
box()

# dev.off()


######## Table 1 is case summary of presidents who also served in Congress 
# it is the data in presModerationData with formatting/editing for better reading


######## Figure 2 plots ideal points as President versus Member of Congress

# Uncomment next line and the corresponding dev.off() call to generate .png file images
# png(file="moderationResultsFigure.png",width=4,height=4, units="in", pointsize=10, res=600)

par(mar=c(3.1,3.1,0.6,0.6), omi=c(0.1,0.1,0.1,0.1), family="serif", mfrow=c(1,1))
plotSize <- .7
plot(x="",y="", xlim=c(-plotSize,plotSize), ylim=c(-plotSize,plotSize), axes=F,xlab="",ylab="")
axis(side=1,at=seq(-plotSize,plotSize,by=.2), line=0, cex.axis=1)
axis(side=2,at=seq(-plotSize,plotSize,by=.2),las=2, line=0, cex.axis=1)
thisLabelY <- expression("Ideal Point as President, president"[i])
thisLabelX <- expression("Ideal Point as Member of Congress, memcong"[i])
mtext(side=1,text=thisLabelX,cex=1,line=2)
mtext(side=2,text=thisLabelY,cex=1,line=2.5)
polygon(x=c(0,2,2), y=c(0,0,2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(0,-2,-2), y=c(0,0,-2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(0,0,2), y=c(0,2,2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
polygon(x=c(0,0,-2), y=c(0,-2,-2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
points(x=dwnom1cong, y=dwnom1pres, pch=pointpch ,col="black",cex=.8)
text(x=dwnom1cong, y=dwnom1pres, label=name, pos=textpos,offset=.3,cex=.6)

legend(x=.20,y=-.5, c("Supports Moderation","Supports Amplification","Not Expected"), 
       col = c("black","black","black"),pch=c(16,15,8),
       text.col = "black",bg = "white",cex = .7,ncol=1, pt.cex=.7, box.lwd=0)
box()
# dev.off()



######## Table 2 reports results of estimating four regression models 

# Model 1
olsModel <- lm(dwnom1pres ~ dwnom1cong)
summary(olsModel)

# Model 2
olsModel2 <- lm(dwnom1pres ~ dwnom1cong + (presIndex >= 32) + dwnom1cong*(presIndex >= 32))
summary(olsModel2)


###################  >>>>> need to locate cong_dwnom1_when_mc data
# 0 isn't the right values, but I'm using it as placeholder for analysis
# the variable is described on page 19, note 20. I located the median ideal point
# in Congress during presidencies, but as of 5/11/2025 I am not sure where I saved
# the median ideal points in Congress while presidents were serving in Congress.
# Sorry, I will update if I find it, but I did this analysis over a decade ago and
# I'm unable to retrace all of my steps now. 
cong_dwnom1_when_pres <- medianDW1Congress
cong_dwnom1_when_mc <- rep(0, nrow(presModerationData)) # placeholder only

# Model 3
olsModel_3 <- lm((dwnom1pres - cong_dwnom1_when_pres) ~ (dwnom1cong - cong_dwnom1_when_mc))
summary(olsModel_3)

# Model 4
olsModel_4 <- lm((dwnom1pres - cong_dwnom1_when_pres) ~ (dwnom1cong - cong_dwnom1_when_mc) + (presIndex >= 24) + (dwnom1cong - cong_dwnom1_when_mc)*(presIndex >= 24))
summary(olsModel_4)


########### Figure 3: side by side comparing the modern and historic era presidents

# Uncomment next line and the corresponding dev.off() call to generate .png file images
# png(file="moderationResultsFigureTimeVarying.png",width=8,height=4, units="in", pointsize=10, res=600)

par(mfrow=c(1,2),mar=c(3.1,3.1,2.1,0.6),omi=c(0.1,0.1,0.1,0.1),family="serif")
plotSize <- .7
plot(x="",y="", xlim=c(-plotSize,plotSize), ylim=c(-plotSize,plotSize), axes=F,xlab="",ylab="",main="(a) Historic Era Presidents")
axis(side=1,at=seq(-plotSize,plotSize,by=.2), line=0, cex.axis=1)
axis(side=2,at=seq(-plotSize,plotSize,by=.2),las=2, line=0, cex.axis=1)
thisLabelY <- expression("Ideal Point as President, president"[i])
thisLabelX <- expression("Ideal Point as Member of Congress, memcong"[i])
mtext(side=1,text=thisLabelX,cex=1,line=2)
mtext(side=2,text=thisLabelY,cex=1,line=2.5)
polygon(x=c(0,2,2), y=c(0,0,2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(0,-2,-2), y=c(0,0,-2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
abline(v=0)

# here is the regression line for historic era only
abline(a=0.017, b=0.424, lty=3, lwd=2, col="black")
points(x=dwnom1cong[presIndex<32], y=dwnom1pres[presIndex<32], pch=pointpch[presIndex<32] ,col="black",cex=.8)
text(x=dwnom1cong[presIndex<32], y=dwnom1pres[presIndex<32], label=name[presIndex<32], pos=c(1,1,3,2,3,1,1,3,1,1,1,1,2,1,4,1,1),offset=.3,cex=.6)

legend(x=-.7,y=.7, c("Supports Moderation","Supports Amplification","Not Expected","Regression Line"), 
       col = c("black","black","black","black"),pch=c(16,15,8,-1),lty=c(-1,-1,-1,3),
       text.col = "black", bg = "white",cex = .7,ncol=1, pt.cex=.7, box.lwd=0)
box()

plot(x="",y="", xlim=c(-plotSize,plotSize), ylim=c(-plotSize,plotSize), axes=F,xlab="",ylab="",main="(b) Modern Presidents")
axis(side=1,at=seq(-plotSize,plotSize,by=.2), line=0, cex.axis=1)
axis(side=2,at=seq(-plotSize,plotSize,by=.2),las=2, line=0, cex.axis=1)
mtext(side=1,text=thisLabelX,cex=1,line=2)
mtext(side=2,text="",cex=1,line=2.5)
polygon(x=c(0,0,2), y=c(0,2,2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
polygon(x=c(0,0,-2), y=c(0,-2,-2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
abline(h=0)

# here is the regression line for modern presidents
abline(a=0.017+0.095   ,b=0.424+0.976   ,lty=3,lwd=2,col="black")
points(x=dwnom1cong[presIndex>=32], y=dwnom1pres[presIndex>=32], pch=pointpch[presIndex>=24] ,col="black",cex=.8)
text(x=c(-0.373,0.217,0.301,-0.370,-0.318,-0.244,0.167),
     y=c(-0.314,0.437,0.452,-0.321,-0.434,-0.370,0.469),
     label=c("OBAMA","BUSH HW","FORD","TRUMAN","KENNEDY","JOHNSON LB","NIXON"),
     pos=c(3,1,3,1,1,1,3),
     offset=.3,cex=.6)
legend(x=.15,y=-.40, c("Supports Moderation","Supports Amplification","Not Expected","Regression Line"), 
       col = c("black","black","black","black"),pch=c(16,15,8,-1),lty=c(-1,-1,-1,3),
       text.col = "black", bg = "white",cex = .7,ncol=1, pt.cex=.7, box.lwd=0)
box()

# dev.off()

# misc. p-value calculation
# dnorm((0.57768 - 1) /    0.11855)

########## Various Robustness Checks

# were historic pres less extreme as members of congress?
mean(abs(na.omit(dwnom1cong[presIndex>=32]))) - mean(abs(na.omit(dwnom1cong[presIndex<32])))
sqrt( (var(abs(na.omit(dwnom1cong[presIndex>=32]))) / length(abs(na.omit(dwnom1cong[presIndex>=32])))) + (var(abs(na.omit(dwnom1cong[presIndex<32]))) / length(abs(na.omit(dwnom1cong[presIndex<32])))) )
mean(abs(na.omit(dwnom1cong[presIndex>=32]))) - mean(abs(na.omit(dwnom1cong[presIndex<32])))

# wilson as first modern
olsModel2 <- lm(dwnom1pres ~ dwnom1cong + (presIndex >= 28) + dwnom1cong*(presIndex >= 28))
summary(olsModel2)

# what if models estimated using absolute values?
olsModel_1A <- lm(abs(dwnom1pres) ~ abs(dwnom1cong))
summary(olsModel_1A)
olsModel_2A <- lm(abs(dwnom1pres) ~ abs(dwnom1cong) + (presIndex >= 24) + abs(dwnom1cong)*(presIndex >= 24))
summary(olsModel_2A)

# what if old estimates of presidents' ideal points used?
olsModel <- lm(olddwnom1 ~ dwnom1cong )
summary(olsModel)
olsModel_2A <- lm(olddwnom1 ~ dwnom1cong + (presIndex >= 32) + dwnom1cong*(presIndex >= 32))
summary(olsModel_2A)

# how extreme would modern governors need to be to change conclusion?
# discussed on page 18
dwnom1pres[49] <- -0.406 # F. Roosevelt () 
dwnom1pres[50] <- -0.491 # Carter
dwnom1pres[51] <- 0.549 # Reagan (), 
dwnom1pres[52] <- 0.564 # G.W. Bush ()  
dwnom1pres[53] <- -0.392 # Clinton
presIndex[49] <- 32 
presIndex[50] <- 39 
presIndex[51] <- 40 
presIndex[52] <- 43
presIndex[53] <- 42
dwnom1cong[49] <- dwnom1pres[49]*1.5
dwnom1cong[50] <- dwnom1pres[50]*1.5
dwnom1cong[51] <- dwnom1pres[51]*1.5
dwnom1cong[52] <- dwnom1pres[52]*1.5
dwnom1cong[53] <- dwnom1pres[53]*1.5
olsModel <- lm(dwnom1pres ~ dwnom1cong )
summary(olsModel)
olsModel_2A <- lm(dwnom1pres ~ dwnom1cong + (presIndex >= 32) + dwnom1cong*(presIndex >= 32))
summary(olsModel_2A)

# here's the Monroe robustness check, discussed on page 19
# james monroe, presIndex==5
olsModelDropMonroe <- lm(dwnom1pres[presIndex!=5] ~ dwnom1cong[presIndex!=5])
summary(olsModelDropMonroe)
olsModel2DropMonroe <- lm(dwnom1pres[presIndex!=5] ~ dwnom1cong[presIndex!=5] + (presIndex[presIndex!=5] >= 32) + dwnom1cong[presIndex!=5]*(presIndex[presIndex!=5] >= 32))
summary(olsModel2DropMonroe)

# can also look at diagonal distance and second dimension
# discussed on pp. 19-20 of the article 
distance2d <- function(x1,y1,x2,y2) {
  WEIGHT <- .4438
  a <- (x1 - x2)
  b <- WEIGHT*(y1 - y2)
  c <- sqrt(a^2 + b^2)
  return(c)
}
diagDistPres <- distance2d(dwnom1pres,dwnom2pres,0,0) 
diagDistCong <- distance2d(dwnom1cong,dwnom2cong,0,0)
diagDistPres <- diagDistPres # - mean(na.omit(diagDistPres))
diagDistCong <- diagDistCong # - mean(na.omit(diagDistCong))
olsModel3 <- lm(diagDistPres ~ diagDistCong)
summary(olsModel3)
olsModel4 <- lm(diagDistPres ~ diagDistCong  + (presIndex >= 24) + diagDistCong*(presIndex >= 24))
summary(olsModel4)


png(file="diagonalDistanceSecondDimension.png",width=7.27,height=4, units="in", pointsize=10, res=600)
par(mfrow=c(1,2),mar=c(3.1,3.1,2.1,0.6),omi=c(0.1,0.1,0.1,0.1),family="serif")
plotSize <- .75
plot(x="",y="",asp=1,xlim=c(0,plotSize),ylim=c(0,plotSize),axes=F,xlab="",ylab="",main="(a) Diagonal Distances")
polygon(x=c(-1,2,2), y=c(-1,2,-1), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(-1,-1,2), y=c(-1,2,2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
points(y=diagDistPres[diagDistPres<diagDistCong],x=diagDistCong[diagDistPres<diagDistCong],pch=16)
points(y=diagDistPres[diagDistPres>=diagDistCong],x=diagDistCong[diagDistPres>=diagDistCong],pch=15)
abline(a=0,b=1,col="black")
mtext(side=1,text="Diagonal Distance as Member of Congress",cex=1,line=2)
mtext(side=2,text="Diagonal Distance as President",cex=1,line=2.5)
axis(side=1,at=seq(0,plotSize,by=.15), line=0, cex.axis=1)
axis(side=2,at=seq(0,plotSize,by=.15),las=2, line=0, cex.axis=1)
box()


plotSize <- 1.0
plot(x="",y="",asp=1,xlim=c(-plotSize,plotSize),ylim=c(-plotSize,plotSize),axes=F,xlab="",ylab="",main="(b) Second Dimension Scores")
polygon(x=c(0,2,2), y=c(0,0,2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(0,-2,-2), y=c(0,0,-2), density = NULL, angle = 45,
        border = NULL, col = "gray90", lty = 1)
polygon(x=c(0,0,2), y=c(0,2,2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
polygon(x=c(0,0,-2), y=c(0,-2,-2), density = NULL, angle = 45,
        border = NULL, col = "gray80", lty = 1)
points(y=dwnom2pres[abs(dwnom2pres)<abs(dwnom2cong) & dwnom2pres > 0 & dwnom2cong > 0],x=dwnom2cong[abs(dwnom2pres)<abs(dwnom2cong) & dwnom2pres > 0 & dwnom2cong > 0],pch=16)
points(y=dwnom2pres[abs(dwnom2pres)<abs(dwnom2cong) & dwnom2pres < 0 & dwnom2cong < 0],x=dwnom2cong[abs(dwnom2pres)<abs(dwnom2cong) & dwnom2pres < 0 & dwnom2cong < 0],pch=16)
points(y=dwnom2pres[abs(dwnom2pres)>=abs(dwnom2cong) & dwnom2pres > 0 & dwnom2cong > 0],x=dwnom2cong[abs(dwnom2pres)>=abs(dwnom2cong) & dwnom2pres > 0 & dwnom2cong > 0],pch=15)
points(y=dwnom2pres[abs(dwnom2pres)>=abs(dwnom2cong) & dwnom2pres < 0 & dwnom2cong < 0],x=dwnom2cong[abs(dwnom2pres)>=abs(dwnom2cong) & dwnom2pres < 0 & dwnom2cong < 0],pch=15)
points(y=dwnom2pres[dwnom2pres > 0 & dwnom2cong < 0],x=dwnom2cong[dwnom2pres > 0 & dwnom2cong < 0],pch=8)
points(y=dwnom2pres[dwnom2pres < 0 & dwnom2cong > 0],x=dwnom2cong[dwnom2pres < 0 & dwnom2cong > 0],pch=8)
abline(a=0,b=1,col="black")
mtext(side=1,text="Second Dimension Score as Member of Congress",cex=1,line=2)
mtext(side=2,text="Second Dimension Score as President",cex=1,line=2.5)
axis(side=1,at=seq(-plotSize,plotSize,by=.2), line=0, cex.axis=1)
axis(side=2,at=seq(-plotSize,plotSize,by=.2),las=2, line=0, cex.axis=1)
box()
legend(x=.20,y=-.65, c("Supports Moderation","Supports Amplification","Not Expected"), 
       col = c("black","black","black"),pch=c(16,15,8),
       text.col = "black",    bg = "white",cex = .7,ncol=1, pt.cex=.7, box.lwd=.7)

dev.off()


